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We express the electromagnetic field propagating in an arbitrary time-independent non-dispersive 
medium in terms of an operator that turns out to be pseudo-Hermitian for Hermitian dielectric and 
magnetic permeability tensors £ and fJ>. We exploit this property to determine the propagating 
field. In particular, we obtain an explicit expression for a planar field in an isotropic medium with 

£ — el and \1 = \l 1 varying along the direction of the propagation. We also study the scattering 
of plane waves due to a localized inhomogeneity. 
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The study of electromagnetic (EM) fields propagating 
in a linear (non-dispersive) medium is one of the oldest 
problems of physics. Most textbook treatments of this 
problem begin with the assumption of harmonic time- 
dependence. In this letter, we present a systematic solu- 
tion of this problem that does not rely on this assump- 
tion and instead makes use of the notion of "pseudo- 
Hermitian operator" that was originally developed in 
[l], 0] to deal with "PT-symmetric Hamiltonians Q . [IH 

Consider the propagation of EM fields in a linear 
source-free medium with time-independent dielectric and 
inverse magnetic permeability tensors £ = £ (x) and 
fi = fJ/'(x), Q. Then the electric and magnetic fields, 
E and B, satisfy Maxwell's equations: 

v-z3 = o, v-s = o, (i) 

B + DE = 0, D-T)H = 0, (2) 

where D := £ E, H := /J, B, a dot means a time- 
derivative, and D denotes the curl operator; DF := VxF 
for any vector field F. Our aim is to solve @ for 
E = E(x, t) and B = B(x, t) in terms of the initial fields 
En ■= E(x,0) and B := B(x,0). We will consider loss- 
less media for which £(x) and fi' (x) are (Hermitian) 
positive-definite matrices for all x G M. 3 . 

We begin our study by expressing Eqs. @ in terms 
of E and B. Taking the time-derivative of the second 
of these equations and using the result in the first, we 

obtain E + U 2 E = where 



Q, 2 := £ D fiD. 
We can easily solve E + fl 2 E = to find 



(3) 



In view of the first equation in the magnetic field is 
given by B(x , t) = B (x) - f* ds DE(x, s), [H[. 

Relation (0| is of limited practical importance, because 
in general its right-hand side involves infinite deriva- 
tive expansions. We can choose the initial fields such 
that they are eliminated by a positive integer power of 
fi 2 . This leads to an infinite class of exact solutions of 
Maxwell's equations with polynomial time-dependence. 
In order to use ^ in dealing with physically more inter- 
esting situations, we will express (JSJ and © as integral 
operators and compute the corresponding integral kernels 
(propagators). [l6| This requires a closer look at ft 2 . 

Let -<F, Gy.= J R3 dx 3 F(x)* -G(x), where F and G are 

vector fields. Then H := {F : M 3 -> C 3 | -< F, F y < oo} 
together with the inner product -<•,•>- form a Hilbcrt 
space. It is easy to check that the curl operator D 
is actually a Hermitian operator acting in Tt. That is 

^F^DG^^^QF, G)~. The same is true about D fi'D. 
But Q 2 is not Hermitian. Its adjoint £! 2 ^, which is de- 
fined by the condition -< F, il 2 G fl'^F, G >-, satis- 

fies n 2 1 = £ n 2 £ . This means that fi 2 : H -> H is an 
£ -pseudo-Hermitian operator Here we view £ as an 

operator acting in H according to (£ F) (x) :— £(x)F(x). 

Indeed, because £ (x) is a positive-definite matrix for all 

x, the operator £ is a positive-definite (metric) operator. 
This in turn implies that it defines a new positive-definite 
inner product that renders f2 2 self-adjoint [l|; letting -< 

F,Gy e :=^F,£Gy,weMd^F, Q 2 Gy e =^n 2 F,Gy e . 
Furthermore, f2 2 may be mapped to a Hermitian opera- 
tor h : H — > H via a similarity transforms [2| . A possible 
choice for h is [6| 



E{x,t) = cos(Clt)E (x) + fT 1 sm(ttt)E (x), (4) 



h := £ Q 2 £ = £ D fID£ 



(7) 



where E Q (x) := E(x, 0) = £ D fi' B (x) and 



cos 



00 (_-\\n 
n=0 V 



n- 1 sin(Ot) := t V — — — >— (t 2 n 2 ) n . 
V ' ^ 2n+l ! v ' 



(5) 
(6) 



Note that because £ (x) and /i (x) are assumed to be 
positive matrices for all x, they have a unique positive 
square root This in turn implies that h is a posi- 
tive operator with a nonnegative spectrum and a unique 
positive square root h?, 171 ]. 

Because h is Hermitian, we can use its spectral reso- 
lution to compute any function $ of h, Q. In light of 



2 



$(f2 2 ) = S 2 £ 2 , this allows for the calculation of 
the action of $(£! 2 ) on any vector field G: 

<f(n 2 )G(x)= [ dy 3 e(x)-i{xmh)\y)E(y^G(y), (8) 

JR 3 

where we have used Dirac's bra-ket notation. 

To demonstrate the effectiveness of the above method 
we consider the textbook Q problem of the planar prop- 
agation of the initial fields 



where 



E (z) = £ (z) e~^ z i, B (z) = B (z)e- l «° z j, (9) 

along the z-axis in an isotropic medium with £ [x] = 

e{z) 1 and jj, (x) = where x =: (x,y,z), 1 

is the identity matrix, £q and Bo are given envelope 
functions, ko is the principal wave number at which the 
Fourier transform of Eq(z) and Bq(z) are picked, i and 
j are the unit vectors along the x- and y-axes, and e(z) 
and fJ-(z) are respectively (the z-dependent) dielectric and 
magnetic permeability constants. We will in particular 
consider the cases that e(z) and tend to constant 
values as z — > ±oo. 

For this configuration all the fields are independent of 
x and y-coordinates and we have Q 2 = e(z)~ l p n{z)~ l p, 
where p := —i-j-, 



e(z) *p(j,(z) 1 pe{z) a. 



(10) 



and Eq{z) := e{z)- l D Q {z) = e{z)~' L 'SB Q {z) = 
eizy^ikoBoiz) - B' (z)] e~^ z l 

In order to determine the spectral resolution of h we 
need to solve the time-independent Schrodinger equation 
for the position-dependent- mass Hamiltonian (|10|) , i.e., 



-eizyi-^Lizr^Uzrhiz)]) =^{z). (11) 

Because of the above-mentioned asymptotic behavior of 
e(z) and /J,(z), the eigenfunctions of h are the solutions of 
pip fulfilling the bounded boundary conditions at ±oo. 
Also note that uj 2 G M + , because h is a positive operator. 

For an arbitrary e we cannot solve (1111) exactly. There- 
fore, we employ the WKB approximation. To do this we 
express ip in its polar representation: ip = R e where 
R and S are real-valued functions. Inserting ip = R e lS 
in (fTTj) gives 



S'(z) 2 + Q(z) =LU 2 e(z)ii(z), 

±Mz)- 1 R-(z) 2 S'(z)]=0, 
dz 



where 

Q{z) 



Ai(z)- 1 i?_(z) ' 



J2_(z) 



e(z)2 



(12) 
(13) 

(14) 



WKB approximation amounts to neglecting Q(z) in (|12|) . 
This yields 

S(z) — lu u(z) + ci, R-(z) = c 2 e(z)~^ [i(z)i , (15) 



u{z):= div{i)-\ v(z):=[e(zMz)]-K (16) 



and ci , C2 are possibly w-dependent integration con- 
stants. Using these choices for S and i?_ and fixing 
ci and C2 appropriately, we find the following (S-function 
normalized eigenfunctions for all u> £ K. 



ipoj(z) 



e iu>u(z) 
^Ttv(z) 



(17) 



Next, we use ip w to express cos(h 2 t) and /i 2 sin(/i 2 i) 
in terms of their spectral resolution: 



cos(/i 2 i) 



cos(a;t) \ip u )(ipu 



(18) 



fe-f sin(/»*f) = I doj hMGM- (19) 



i . . , i 



Using (H7D - dUD and the identities J^dcu e lau3 = 
2tt5(o) and /_ dcu e mw /u — iri sign(a), with 5 denoting 
the Dirac delta function and sign(a;) := x/\x\ for x ^ 
and sign(0) := 0, we find 



(z\ cos(h 2 t)\w) 



dw cos(wi) ipu(z) ip u (w)* 

-oo 

[v(z)v(w)]~*A(z,w;t), (20) 



(z\h 2 sm(h 2 t)\w) = / c&ja> 1 sin(a;i) ip u (z) ip u (wy 



where 



= -[v(z)v(w)]-^(z,w;t), (21) 



A(z, t) :— S[u(w) — u{z) + 1] + S[u(w) — u(z) — t], 

S(z, w; t) :— sign[w(u>) — u{z) +t] — sign[u(w) — u(z) — t\. 

Because u is a monotonically increasing function that 
vanishes only at z = 0, it is invertible and its inverse u 
is also a monotonically increasing function with a single 
zero at z = 0. This implies that the quantity 

w±(z,i) := u _1 (m(z) ±<) (22) 

is the only zero of u(w) — (u(z) ± t). Hence, 

r r / \ / s , .1 S[W - W T (z,t) 

6[u( W )-u(z)±t]= WMzM - 
In view of this relation and (1161). we have 



a / ,n , S[w-w + (z,t)] 

A(Z,wt) = ; ; rr : ; ; — — . (2d) 

v ' ' ' v{w^{z,t))- 1 v{w + {z,t))- 1 v ; 
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Furthermore, because both u and u 1 are monotonically 
increasing, for t > we have W-(z,t) < w+(z,t) and 



E(z, w; t) 



for W-(z,t) < w < w + (z,t), 

otherwise. ^ ' 



Next, we compute the action of cos(rii) on an arbitrary 
test function f(z). To do this we set $(f2 2 ) = cos(f2i) in 
(|8j and use ([20]) and (p3|) to evaluate the corresponding 
integral. This yields 



CQB(nt)f(z) = f-(z,t) + f+(z,t), 



(25) 



where 



f±(z,t) := 



1 



e(w±(z, t))v{w±(z, t)) 



Similarly, setting $(S! 2 
(dU and dMI) we find 



e(z)u(z) 
) = fi- 1 sin(Oi) in 



dSJ and using 



fT 1 sin(fit)/(z) = 



/ 






v(w) 5 



Finally, we use tfTS]) . and (gSJ) to express (g]) as 



(26) 



£(z,t) 




e(to_(z,t)) 



/i(u)+(z,i))_ 
+(*,*) 



1 
4 



£*)(«>_(*,*)) + 



-(*>t) 



dw; [i(w) 4 e(w) 4 Eq(w) 



(27) 



According to (|22|) . u>±(z,0) = z. This shows that 
setting t = 0, the integral on the right-hand side of 
(|27p disappears and the remaining terms add up to give 
Eq(z). Similarly we can use (fT6|) and (f22|) to establish 
^±(Z)0) = ±v(z) = ±[e(z)/Li(z)] _ 2 and use the latter 
to show that setting t = in the time-derivative of the 

right-hand side of (j27|) yields E a {z). 

In vacuum where v = (e/i) -1 ' 2 = c, WKB approxi- 
mation is exact, u(z) — z/c, w±(z,t) — z ± ct, and (127|) 
gives 



i{£ (z-ci)+£ (2 + rf)+ 

fZ-\-Ct . ^ 

du; £„(») I , (28) 



£(z,i) 



which is precisely D'Alembert's solution of the 1+1 di- 
mensional wave equation [ic| |. 

In general, (|2"T|) is a valid solution of Maxwell's equa- 
tions, if WKB approximation is reliable. This is the case 
whenever |Q(z)j is negligibly smaller than the right-hand 
side of (fT2"]) . [18l |. In view of (fT5|) and (jTHJ) this condition 
takes the form 



v 

Y 



vv 



1 /2 
— V 

2 _ 



< uj 2 , 



(29) 



where we have suppressed the z dependence of v and fx. 
Due to the asymptotic behavior of \x and e, v tends to 
constant values as z — > ±oo. This in turn implies that the 
square root of the left-hand side of (f2"9"]) has a least upper 
bound that we denote by u> m in- In this case, (|29[) means 
\uj\ 3> Wmin- Recalling the role of to in our derivation 
of (|27[) . we can view this condition as a restriction on 
the choice of the initial conditions. More specifically, 
(127[) is a good approximation provided that for all u> G 

[-Wmin, Wmin], (V'w |e* |#0 (z) |) ~ and (^w|e'|-Eo(z)|) ~ 
0. For a planar laser pulse with initial envelope functions 
So and Bq picked far away from the region where s and 
li vary significantly, these conditions hold for c u),^ <C 
|fco|. The same is true for an initial plane wave with 
sufficiently large wave number \ko\. 

Next, we wish to elaborate on the relationship between 
our approach and the standard application of the WKB 
approximation in solving Maxwell equations, particularly 
for the effectively one-dimensional model we have been 
considering. The first step in the direct application of the 
WKB approximation to Maxwell equations is to assume 
a harmonic time-dependence for the solution. This yields 
an ordinary differential equation that can be treated us- 
ing the WKB approximation, Q. If we apply e -1 / 2 to 
eigenfunctions (|17[) of the Hermitian operator h we ob- 
tain (WKB-approximate) eigenfunctions of the pseudo- 
Hermitian operator Q 2 , Q. It is very easy to check that 
these are solutions of the Maxwell equations ([2]) with a 
harmonic time-dependence. Up to a trivial normaliza- 
tion constant, they coincide with the conventional WKB 
solutions: 



e(z) 



Furthermore, if we make the additional assumption that 
form a complete set in the sense that for any solution 
E(z,t) of Maxwell's equations, there are c(u>) 6 C such 
that 



E(z,t) = / du c{uj)E w {z,t), 



(30) 



then we can use the initial electric and magnetic field to 
determine c(u>) and use the result to evaluate the integral 
in ([30]) . The end result of this lengthy calculation is 
identical to ([27)1 . 

The observations that fl 2 is pseudo-Hermitian and that 
it can be mapped by a similarity transformation to a Her- 
mitian operator h provide a simple justification for the 
above-mentioned assumptions regarding harmonic time- 
dependence and completeness of the conventional WKB 
solutions. The consequences of the pseudo-Hermiticity 
of £1 2 is not limited to the application of the WKB ap- 
proximation. For instance, it implies the existence of ex- 
act solutions of Maxwell's equations that have harmonic 
time-dependence and form a complete set. 

As a concrete example of the application of (12~7| , sup- 
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pose that s has a Lorentzian shape and p, is a constant: 



s(z) =e Q l + a(l + 7 2 z 2 ) 1 



M = Mo, (31) 



where e , a, 7, /i n are positive constants and £$Pq = c 
Then, by inspection, we can show that 



< 



c^3a(l + va) 
27 ' 



where v := 1 + 45/256 « 1.176. For a < 1, w min < C7 _1 . 
This means that WKB approximation and consequently 
(|2T|) are valid provided that the allowed w values be much 
larger than C7 _1 . 

In order to implement (|27|) , we must compute w± . This 
involves the evaluation and inversion of u. For the choice 
(|3ip . we expand u and u" 1 in power series in the inho- 
mogeneity parameter a and perform a perturbative cal- 
culation of w±. This yields 



w±(z,t) = z±ct=fa 



70±(z,t) 



± 



2/ 7 
1 16 



{\(z,t)6 ± (z,t) + v ± (z,t)}}+0(a 3 ), (32) 



where 



0±(z,t) := tan 1 

A(z,t) := 1 + 
v±(z,t) 



7ci 



7 2 + z{z ± ci) 
4 7 2 



(Z ± ct) 2 ' 



7ct[7 2 — z(z ± ct)] 



(7 2 + z 2 )[7 2 + (z± ct) 2 

In view of (f3"Tj) and (f3"2")l , we also have 
1 

e(w±(z, t)) 

= J- -t- « 

3 



= lTtt^,t)±a 2 |((z,f) + 

C(M)T^(M) 2 }+0(a 3 ), (33) 



where 



7 2 ct(2z ± c<) 



4( 7 2 + z 2 )[ 7 2 + (z±ct) 2 ] 

7 3 (z±ct) 
4[ 7 2 + (z±ct) 2 ] 2 ' 



£(*,*) 



With the help of (J32D and Q we can use (21]) to deter- 
mine the dynamical behavior of the EM fields for initial 
configurations of the form ^ provided that the initial 
fields do not violate the condition of the reliability of the 
WKB approximation and that we can neglect the third 
and higher order contributions in powers of a. For ex- 
ample we can use this method to determine the effect of 
the inhomogeneity (|31l) on the planar propagation of a 
Gaussian laser pulse [19j |. 
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FIG. 1: Plots of Ar (full curve) and Ap (dotted curve) as a 
function of z for a = 0.2 in units where 7=1. 



Another application of our results is in the solution of 
the scattering problem. It is not difficult to see that for 
t — > 00 (i.e., ct 3> 2,7), 



w±(z,t) - 
e(w±(z,t)) 



s(z) 



z±ct + Ar(z) for 
±ct ± Ar(0) for z = 

1 

- Ap(z), 



(34) 



(35) 



where 
Ar(z) 

Ap(z) 



a 7, -1 n 

tan 1 

2 

2 

16 

£0 



16 L \z) 



ryz 



7 2 + z 2 



e(z) 
1 - 



07 



+ 0(a 3 ) 
2 5a 2 7 4 



+ 



4( 7 2 + z 2 ) 32( 7 2 + z 2 )< 



+ 0(a 3 ), 



+ 0(a 3 ). 



According to ([271) , (|34|) and ([35| , the scattering of an 
initial plane wave, with Eo = e~ lk ° z i and Bo = 0, by the 
inhomogeneity pip results in a change in the amplitude 
and phase angle of the wave that are respectively given 
by Ap(z) and — k Ar(z). Specifically, as t — > 00, 

E(z, t) -» EJz, t) := E{z, t) Ap(z) e -*> Ar ( z ), (36) 

where E(z,t)\ a=Q := \ (e^M*-*) + e -*k (z+ct)y L The 
predictions of (|36p should be experimentally verifiable in 
typical interferometry experiments [20| . Figure [1] shows 
the plots of Ap and Ar for a = 0.2. As seen from this 
figure Ar has a discontinuity at z = 0. 

A quantity of direct physical relevance is the Fourier 
transform E s (k,t) := ± J^dz e lkz E s {z,t) of E s (z,t). 
Up to linear terms in a, it is given by 

E s (k, t) = \s(k -k) + aAE(k a ,k) + 0(a 2 )] cos(fc ci)i, 



where AE(ko, k) := 7 



l-e 



-~,\k -k\ 



-~f\k -k\ 



Figure [2] 



4(l-fc/fc ) 8 

shows the graph of AE(ko, k) as a function of k. Again 
there is a discontinuity at k = fep. 
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k 




0.998 1 1.002 

FIG. 2: Plot of AE(k ,k) as a function of k for 7 = 10" 3 m 
and ko — 10 7 m _1 . fc is measured in units of fco = 10 7 m _1 . 

To summarize, we have obtained a closed form expres- 
sion ([4| for the propagation of EM fields in an arbitrary 
non-dispersive stationary medium that yields the fields in 
terms of a pseudo-Hermitician operator fl 2 . This allows 
for formulating the problem in terms of an equivalent 
Hermitian operator h. Using the spectral resolution of 
h and the WKB approximation we have found an ex- 
plicit formula, namely 1(27]). for the propagating fields 
and demonstrated its application for the scattering of 
plane waves moving in an inhomogeneous non-dispersive 
medium. Although similar spectral techniques have pre- 
viously been used in dealing with EM waves [TTI | . we 



believe that our treatment provides a more straightfor- 
ward and systematic solution for this problem. Indeed 
to our knowledge an expression for the propagating field 
as general and explicit as Eq. (|27[) has not appeared in 
the literature previously. Another important aspect of 
our method is the wide range of its applications, e.g., it 
can be used to study wave propagation in inhomogeneous 
fibers, observation of superfluid vortices, etc. 

Our results may be generalized in various directions. 
For example, for the cases that £ fails to be Hermitian, 
one may appeal to the notion of weak pseudo-Hermiticity 
[T3 | and use the results of [HI to obtain an appropriate 
equivalent Hermitian operator h to fl 2 . One may also 

incorporate the dispersion effects by letting the £ and 
fi that appear in the eigenvalue equation for h to de- 
pend (via a dispersion relation) on u>. This will lead to a 
modification of (|27p that we plan to explore in future. 
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